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Abstract: Estimates for the sedimentation rate of realistic ice crystals at sizes smaller than 100 /im are presented. These calculations, 
which exploit new results for the capacitance of ice crystals, are compared with laboratory studies and found to be in good agreement. 
The results highlight a weakness in contemporary ice particle fall speed parameterisations for very small crystals, which can lead to 
sedimentation rates being overestimated by a factor of two. The theoretical approach applied here may also be useful for calculating 
the sedimentation rate and mobility of non-spherical aerosol particles. 
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1 Introduction 

The sedimentation of natural ice particles is fundamen- 
tal to the evolution of the ice phase in clouds, influencing 
almost all of the relevant microphysical processes includ- 
ing deposition, riming, aggregation, evaporation and melt- 
ing. Modelling the sedimentation velocity v of ice crystals 
is challenging because they are almost universally non- 
spherical in shape, and span a range of flow regimes. 
Whereas a rather accurate formula exists for the drag on 
a sphere (Abraham 1970), theoretical progress for non- 
spherical particles is more difficult. Importantly, Bohm 
(1989) and Mitchell (1996) have attempted to modify 
Abraham's formulation for an ice particle of arbitrary 
shape: Mitchell argues that the results are likely to be 
accurate to within ~ 20% over the whole range of flow 
regimes experienced by ice particles in the atmosphere. 
His piecewise fits have recently been refined to pro- 
vide a continuous power law with variable co-efficients 
(Khvorostyanov and Curry 2002), and further adjustments 
have been made for turbulence at high Reynolds number 
(Mitchell and Heymsfield 2005, Khvorostyanov and Curry 
2005). This method has proved rather successful, provid- 
ing fall velocities consistent with observations for a range 
of size and shape particles, and the framework itself is 
very appealing because of its generality. However there is 
a need to test it, particularly for very small crystals where 
the boundary-layer ideology underpinning it breaks down. 

Ice crystals smaller than 100/zm can play an impor- 
tant role in clouds. Such crystal sizes are common at the 
top of stratiform clouds and Rauber and Tokay (1991) 
suggest that it is this feature which allows supercooled liq- 
uid droplets to persist there in a weak updraught, despite 
the flux of vapour to the ice crystals. The sedimentation 
of the crystals is key to a quantitative understanding of 
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this situation, since the longer the crystals reside in the 
mixed-phase layer, the more time they spend growing at 
the droplets' expense. This range of crystal sizes is also 
typical of supercooled fogs (eg. Yagi 1970) and diamond 
dust. Aircraft contrails are composed of crystals of order 
lO^rn in size, and their sedimentation influences the extent 
to which the contrail removes moisture from the upper 
troposphere, and could affect its persistence/transition to 
cirrus. In natural cirrus, there remains much controversy 
over the concentrations of ice crystals smaller than around 
60^m. Aircraft measurements have indicated that these 
tiny crystals may exist in very large numbers (see Heyms- 
field and McFarquhar 2002) and could dominate the opti- 
cal properties of the cloud; on the other hand, there is also 
evidence to suggest that some of the measurements of this 
small crystal mode are an artefact of larger crystals shat- 
tering on the probe inlet (Field et al. 2003, McFarquhar 
et al. 2007, Heymsfield 2007). Doppler lidar observations 
could inform this debate, since if these tiny crystals were 
to dominate the optical properties of the cloud, they should 
also dominate the lidar backscatter and Doppler veloc- 
ity. Given accurate estimates of the sedimentation rates of 
these small ice crystals, there is an opportunity to compare 
these values with the measured frequency distribution of 
lidar Doppler velocities in cirrus clouds, and to determine 
whether or not sub-60^m crystals genuinely have a signif- 
icant impact on the optical properties of cirrus. 

The aim of this short paper is to investigate the sed- 
imentation velocity of small ice crystals where the air 
flow is dominated by viscous forces. New estimates of 
the 'capacitance' for realistic ice crystal shapes (West- 
brook et al. 2008) allows the application of two theoretical 
results from the physics literature (Roscoe 1949, Hub- 
bard and Douglas 1993) to calculate crystal fall speeds; 
these estimates are then compared to experimental data, 
and to the fall-speed parameterisations of Mitchell (1996), 
Mitchell and Heymsfield (2005) and Khvorstyanov and 
Curry (2002, 2005); the latter four studies are collectively 
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referred to as MHKC from now on. Upper and lower 
bounds are also constructed where possible using analyt- 
ical results for spheres and spheroids, providing further 
insight. 

2 Viscous sedimentation speeds 

When the Reynolds number Re = vD/uk is small, viscous 
forces dominate the flow. Here D is taken to be the max- 
imum dimension of the crystal, and vu is the kinematic 
viscosity of air. The Stokes solution for a sphere in a vis- 
cous flow is well known: 



\6nrjJ R 



(1) 



where m and R are the particle mass and radius respec- 
tively, g is the acceleration due to gravity, and 77 is the 
viscosity of the air. Dimensional analysis suggests that an 
identical formula should exist for non-spherical ice par- 
ticles falling in a given orientation (or ensemble of ori- 
entations), but with an 'effective' hydrodynamic radius in 
place of the sphere radius. The problem then is to deter- 
mine this hydrodynamic radius for natural ice crystals. 

Very small crystals are affected by rotational Brown- 
ian motion and often have an approximately random ori- 
entation in the absence of an electric field (eg. Foster and 
Hallett 2002). Hubbard and Douglas (1993) have provided 
a theoretical treatment for the average viscous drag on a 
particle which experiences all possible orientations. They 
preaverage the Oseen tensor (which acts as a Green's func- 
tion for the Navier-Stokes equations) over all angles, and 
find that it is identical to the Green's function for Laplace's 
equation (to within a constant factor). Based on this obser- 
vation they find that the angle-averaged flow system can 
be described through a scalar potential which is fixed on 
the surface of the settling particle, and zero far from it. 
This situation is directly analogous to diffusion and elec- 
trostatic problems. An overview of their analysis is given 
in appendix A, but the end result of this approximation is 
that the hydrodynamic radius for the particle is: 



R = C 



(2) 



where C is the capacitance of the particle in length 
units. The capacitance is commonly used to estimate the 
growth/evaporation of ice particles when the mass flux 
is limited by diffusion of vapour onto the ice surface, 
and it has recently become possible to calculate this par- 
ameter accurately using a Monte Carlo method (West- 
brook et al. 2008). The link between molecular diffusion 
and viscous drag seems sensible since the drag in vis- 
cous flow is essentially the product of momentum dif- 
fusing away from the falling crystal. Inserting the capaci- 
tance for a spheroid into equation 2 yields the exact Perrin 
(1934, 1936) formula; comparison with experimental data 
on non-spherical shapes has yielded agreement to within a 
few percent (Hubbard and Douglas 1993). Blawzdziewicz 
et al. (2005) has provided further theoretical arguments 
to explain why this approximation turns out to be rather 



accurate. Although there is an increasing torque to orient 
the ice crystal as the Reynolds number becomes larger, we 
might expect that the randomly oriented approximation 
may still be acceptable in a number of cases, especially 
for particles which do not have strong symmetry (eg. poly- 
crystals) or column/needle crystals where the orienting 
torque is relatively weak (Katz 1998). 

The strongest orientational torque is likely to be for 
plate-like particles which are observed to orient approx- 
imately horizontal if the crystal grows large enough. For 
example Katz (1998) estimated that a disc with a diameter 
of 60^m and a thickness of 3^m would orient approxi- 
mately horizontally with an angular dispersion of 0.6°. 
In this situation the angle-averaged Hubbard-Douglas 
approximation is likely to underestimate the drag. How- 
ever, Roscoe (1949) provides an approximation for the 
drag on a horizontally oriented particle in viscous flow. 
For a thin planar crystal, the Navier-Stokes equations can 
be transformed into Laplace's equation with a fixed poten- 
tial at the particle surface. Again, this is directly analagous 
to the vapour diffusion problem and the hydrodynamic 
radius is again simply related to the capacitance: 



R 



-C. 



(3) 



The derivation is straightforward, and is reproduced in 
appendix B. Strictly this result applies only to planar 
particles of zero thickness; however comparison between 
the exact formula for a thin oblate spheroid with aspect 
ratio 0.1 (Happel and Brenner 1965), and equation 3 with 
the associated capacitance (C — 0.338Z3 Pruppacher and 
Klett 1997, page 547) reveals a difference of less than 
0.5%. 

The MHKC parameterisations all collapse to the 
same formula in the small Reynolds number limit, and 
predict a hydrodynamic radius proportional to the ratio 
of the projected area of the particle A to its maximum 
dimension D: 



R 
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where Co, Sq are dimensionless parameters characterising 
the boundary layer shape and thickness (which have no 
physical meaning in a viscous flow). For liquid drops 
the parameters suggested by Abraham are used (Co = 
0.292, So = 9.06) which matches the Stokes formula for 
a spherical particle (substituting A = nR 2 , D — 2R). For 
ice particles, MHKC prefer Co = 0.6, £0 = 5.83, giving 
a fall velocity ^15% higher than the Stokes result for a 
spherical particle^. Note that (4) may also be expressed as 
R = (Cq5 2 /24) x 7-D/2 where 7 is the 'area ratio' of the 
particle (ie. the projected area divided by the area of an 
enclosing circle). 

t Khvorostyanov and Curry used the smooth sphere parameters for ice 
particles in their 2002 study, but their subsequent paper (Khvorostyanov 
and Curry 2005) uses C = 0.6, S = 5.83, as do Mitchell (1996) and 
Mitchell and Heymsfield (2005). In this paper we will assume these 
latter parameters. 
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2. 1 Planar crystals 

Planar crystals are a frequent occurence in the atmosphere. 
The standard crystal habit diagram suggests these crystals 
form at temperatures between approximately to — 3°C 
and — 8°C to — 25° C; however recent laboratory exper- 
iments by Bailey and Hallett (2002, 2004), along with 
in-situ evidence cited in their papers, indicates that this 
diagram is probably biased by the use of silver iodide 
in laboratory studies, and plate-like crystals can form at 
temperatures as cold as — 70°C, particularly at low super- 
saturations. 

2.1.1 Hexagonal plates 

The capacitance of a hexagonal plate is given by West- 
brook etal. (2008): 

C = 0.58a x [1 + 0. 95(i/2a) 75 ] (5) 

for a crystal of span (2a) across the basal face, and 
thickness L, with aspect ratio L/2a < 1. The fall speeds 
of crystals with aspect ratios in the range 0-0.7 have 
been calculated using the Hubbard-Douglas and Roscoe 
approximations: these curves are shown in figure 1. To 
facilitate the comparison between theory and experiment 
the fall speeds have been normalised relative to the fall 
speed of a sphere with the same mass and maximum 
dimension (ie. normalised velocity^ D/2R). 

Experimental data comes from Michaeli (1977) and 
Kajikawa (1973). Michaeli measured the mass, dimen- 
sions, and sedimentation velocities of free-falling hexag- 
onal plate crystals in a cloud chamber with D ranging 
from 50 to lOO^m. From this information R and the nor- 
malised velocity were calculated directly: see triangle data 
points in figure 1 . The experimental data lies between the 
Roscoe and Hubbard-Douglas curves, somewhat closer to 
the former, suggesting the crystals were to a large extent 
oriented in the horizontal. The point marked by an asterisk 
appears to be an outlier; note that for this particular point 
the thickness of the crystal was not directly measured, but 
calculated using an empirical relationship from a different 
study - the measured velocity would be more consist with 
an aspect ratio closer to unity. Similarly, Kajikawa (1973) 
seeded a supercooled fog in a cold room and measured the 
fall speeds of the resulting crystals (15 — 100/im)*. Fall 
streak photographs were used to calculate the velocities, 
and mass and thickness measurements were taken from 
similar crystals of the same type and diameter which fell 
out at the same time. There is a large amount of scatter in 
his data; however, taking the bin averages from his table 1 
and comparing them with his calculated values for circular 
discs the associated normalised velocties can be estimated 
(diamonds). The aspect ratios were estimated from his fig- 
ure 5. Again the data points lie between the Roscoe and 

* Kajikawa's 20/jm bin was not used here because the fall speeds were 
measured to be higher than that of an inscribed disc falling edge-on, 
suggesting that the melted diameter of these tiny crystal is likely have 
been underestimated. 
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Figure 1. Fall speed of hexagonal plate ice crystals, normalised 
by the fall speed of a sphere with the same mass and maximum 
dimension. Solid black line shows Hubbard & Douglas prediction 
for randomly oriented crystals (R — C), dashed black line shows 
the Roscoe prediction for horizontally oriented crystals (R = | C). 
Diamonds and triangles show fall speeds of ice crystals as measured 
by Kajikawa (1973) and Michaeli (1977) respectively. Asterisk indi- 
cates the aspect ratio for this crystal was not measured (see text). 
Hollow circles are measurements of circular discs falling horizon- 
tally. Filled circles at L/2a — represent analytical predictions 
for discs of zero thickness falling horizontally (black) and edge-on 
(grey). MHKC predictions are shown by the grey solid line (ran- 
domly oriented) and grey dashed line (horizontally oriented). Shad- 
ing indicates fall velocities which lie outside the upper and lower 
bounds derived from analytical results (see text). 

Hubbard-Douglas curves indicating horizontal orientation 
with a certain amount of flutter. 

Jayaweera and Ryan (1972) also measured the fall 
speed of small plate-like crystals; however, unlike their 
measurements of columnar crystals (see section 2.2) they 
did not explicitly record the crystal aspect ratio. Their fit 
corresponds to a normalised fall speed of 1.41; assuming 
an aspect ratio of ~ 0.1 this value lies between the Roscoe 
and Hubbard-Douglas curves. 

Also shown in figure 1 are analytical results for 
infinitely thin circular discs, and experimental data for 
discs of finite thickness (see Clift et al. 1978). These discs 
have the same diameter D as the hexagonal plates. For a 
thin disc falling horizontally the drag is slightly larger than 
for the hexagonal plate prediction, with the disc falling 
9% slower than the plate. For larger aspect ratios the hor- 
izontal discs fall slightly faster than the plates, indicating 
that the Roscoe curve underestimates v somewhat as we 
depart from the zero thickness approximation. For a thin 
circular disc falling edge-on the normalised velocity is 
slightly higher than the Hubbard-Douglas prediction for a 
randomly oriented hexagonal plate, which seems sensible. 

The MHKC formula was used to calculate the grey 
curves in figure 1 assuming both randomly oriented 
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(A = a 2 [1.3 + 3(L/2a)]) and horizontally oriented crys- 
tals (A — 2.6a 2 ). The former assumption is recommended 
by Mitchell (1996) for crystals with Re < 10. The curve 
for randomly oriented crystals predicts fall speeds approx- 
imately 70-80% higher than the experimental data for 
aspect ratios between 0.1-0.3 where most of the data 
points are clustered. The curve for horizontal orientation 
gives more realistic fall speeds for thin plates, around 10- 
20% higher than the experimental data. The Jayaweera 
and Ryan data fit is in approximate agreement with the 
horizontal prediction if the aspect ratio of the crystals 
is assumed to be around 0.1. The strong sensitivity of 
eqution 4 to crystal orientation appears to be unphysi- 
cal given the relatively weak dependence on fall attitude 
observed for circular discs (for very thin discs edge-on is 
45% faster than for horizontal orientation, a dependence 
that becomes weaker for increasing aspect ratio, see Clift 
etal. 1978). 

Hill and Power (1956) showed that bounds could be 
constructed for the drag on a particle of arbitrary shape 
by using known results for inscribed and circumscribed 
spheroids. In this vein the fall speed for an oblate spheroid 
inscribing the hexagonal plate crystal was calculated using 
Oberbeck's formula (Happel and Brenner 1965). The 
edge-on orientation was identified as that corresponding 
to the maximum possible fall velocity. Any values of v 
larger than this upper bound must therefore be erroneous, 
and this indicated by grey shading in figure 1 . The MHKC 
curve for randomly oriented crystals violates this bound 
for aspect ratios less than 0.3; similarly the horizontally 
oriented curve violates the bound for aspect ratios larger 
than 0.55. The Hubbard-Douglas and Roscoe curves both 
lie below the upper bound for all aspect ratios. A second 
bound may also be constructed (not shown in figure 1 for 
clarity) for plates which are assumed to lie in a horizontal 
orientation. The MHKC curve for horizontal plates lies 
above this upper bound by 10-20% for all aspect ratios; 
the Roscoe curve lies below it for all aspect ratios. 

It is also possible to construct a lower bound. 
Using the drag for a circumscribed sphere (diameter = 
\/i 2 + (2a) 2 ) such a bound for was constructed: value 
lower than this bound are shaded in figure 1 . The Roscoe 
curve dips slightly below this lower bound for aspect 
ratios larger than 0.25, as do a couple of the experimen- 
tal data points (by a couple of percent). For an aspect 
ratio of 0.7 the Roscoe curve is 8% lower than the enclos- 
ing sphere lower bound. The suggestion is that the zero- 
thickness Roscoe approximation progressively underes- 
timates the fall speeds for increasing aspect ratios. The 
Hubbard-Douglas and MHKC curves lie well above the 
lower bound. 

As the crystals grow larger they may start to devi- 
ate slightly from the viscous drag regime, and this could 
potentially affect some of the experimental data. Keller 
et al. (1967) have shown that in this case the viscous 
drag calculations underestimate the true drag, and there- 
fore the corresponding normalised velocity is reduced. 
Experimental and numerical results for spheres and discs 
suggest that for Re < 1 this deviation is less than 10% 
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(Clift et al. 1978). The experimental data for ice crystals 
all had Reynolds numbers Re < 0.4. 

Note that there is a slight subtlety in the meaning 
of the maximum dimension D in the above analysis: 
for planar crystals we take D to mean the maximum 
dimension across the basal face ie. D = 2a since this is 
the dimension which is usually measured in observations. 
If D is taken to mean the three dimensional maximum 
span (= y/L 2 + (2a) 2 ) then all the normalised velocities 
are shifted up by a factor / = y/l + (L/2a) 2 except for 
the MHKC curves which are shifted up by a factor f 2 
since D also appears in the definition of R via equation 
4. The result of this transformation is that the entire 
MHKC curve for randomly oriented crystals now lies 
above the upper bound derived above, as do the MHKC 
predictions for horizontally oriented plates with L/2a > 
0.25; the differences with the experimental data points are 
magnified by a factor /. 

2.7.2 Branched & dendritic crystals 

Kajikawa (1973) and Michaeli (1977) also made measure- 
ments of the fall speed of branched and dendritic crystals, 
and their experimental results are shown in figure 2. Based 
on photographs in Kajikawa's paper the capacitance of a 
'representative' branched crystal has been calculated for 
aspect ratios between L/2a = 0.05-0.35 using the method 
outlined in Westbrook et al. (2008). The outline of the 
model crystal is shown inset in figure 2. The width of 
the branches is 35% of their length from the crystal cen- 
tre. The overall span from tip to tip across the crystal is 
defined as (2a) and the thickness as L in the same way 
as for the hexagonal plate crystals. The values of C were 
only marginally lower than for a solid hexagonal plate. 
The capacitance data were substituted into equations 3 and 
2, and curves fitted to the resulting fall speed predictions: 
these are shown in figure 2. The maximum dimension was 
defined as the span from tip to tip across the projection 
shown in figure 2, ie. D = 2a. 

The experimental data appear to be in broad agree- 
ment with the range of velocities indicated by the Roscoe 
and Hubbard-Douglas theories, although the scatter is 
rather wider than for the hexagonal plate data, proba- 
bly indicating the range of crystal shapes and perhaps 
a broader distribution of crystal orientations. The aspect 
ratio of the crystals observed by Michaeli (triangles) were 
not directly measured but estimated using an empirical 
relationship from a different study, further contributing to 
the experimental scatter. 

Also shown in figure 2 are the MHKC predictions for 
the model crystal described above. Under the horizontal 
orientation assumption the MHKC curve is around 45% 
higher than the experimental data, whilst for the randomly 
oriented assumption the MHKC curve is around 80% 
higher. This is particularly an issue for crystals with 
thin branches and dendritic features where the projected 
area A becomes very small but the capacitance is only 
marginally reduced relative to a solid hexagonal plate 
(typically by only ~25% even for quite thin branches, 
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Figure 2. Normalised fall speed of branched and dendritic planar 
ice crystals. Solid black line shows Hubbard-Douglas prediction 
for randomly oriented crystals (R — C), black dashed line shows 
the Roscoe prediction for horizontally oriented crystals (R = |C). 
Diamonds and triangles show fall speeds of ice crystals as measured 
by Kajikawa (1973) and Michaeli (1977) respectively. MHKC 
predictions are shown by the grey solid (randomly oriented) and 
grey dashed line (horizontally oriented). Inset is the branched crystal 
model used in the capacitance and projected area calculations. 
Shading indicates fall velocities which lie outside the upper and 
lower bounds derived from analytical results (see text). 



see Westbrook et al. 2008). If the capacitance is the 
appropriate length scale then this implies that MHKC 
theory will increasingly overestimate the crystal fall speed 
as the crystals become more tenuous. 

It was not possible to construct an upper bound 
using an inscribed oblate spheroid as before, because of 
the branched nature of the crystal. However, Weinberger 
(1972) showed that a sphere with equal volume acts as an 
upper bound for the settling velocity of a non-spherical 
particle in a viscous flow. The normalised velocity for a 
sphere with the same volume as the model crystal was cal- 
culated as shown by the dark grey region in figure 2. The 
MHKC curve for random orientation lies entirely above 
this upper bound. The curve for horizontal orientation 
also violates the upper bound for aspect ratios larger than 
1/8. The Roscoe and Hubbard-Douglas approximations lie 
below the upper bound for all aspect ratios. 

A lower bound may be constructed in the same way 
as for plate crystals and this is indicated by the shaded 
region in figure 2. A couple of the experimental data points 
from Michaeli fall below this bound, although given the 
uncertainty in the aspect ratio this could be an artefact for 
the crystal plotted at L/2a = 0.25. All of the theoretical 
curves lie above the lower bound. 
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Figure 3. Normalised fall speed of hexagonal column and needle 
ice crystals. Solid black line shows Hubbard-Douglas prediction for 
randomly oriented crystals (R — C). Stars, diamonds and triangles 
show fall speeds of ice crystals as measured by Jayaweera and Ryan 
(1972), Kajikawa (1973) and Michaeli (1977) respectively. Circles 
are measurements of circular cylinders falling horizontally (black 
circles) and end-on (grey circles). MHKC predictions are shown 
by the grey solid line (randomly oriented) and grey dashed line 
(horizontally oriented). Shading indicates fall velocities which lie 
outside the upper and lower bounds derived from analytical results 
(see text). 



2.2 Columnar ice crystals 

Hexagonal column and needle crystals typically grow 
at temperatures colder than — 22° C and also in a win- 
dow between -3 and -8°C (Pruppacher and Klett 1997). 
Jayaweera and Ryan (1972) made direct measurements of 
crystal mass, length, width and fall speed; similar mea- 
surements were made by Michaeli (1977), and Kajikawa 
(1973). These data are plotted in figure 3, and include solid 
and hollow hexagonal columns, as well as needles. 

The capacitance of a hexagonal column of length L 
and width (2a) is given by equation 5 with L/2a > 1. Sub- 
stituting this into the Hubbard and Douglas formula (2) 
yields the black line plotted in figure 3. For columns we 
take the maximum dimension to be D = y/L 2 + (2a) 2 . 
Although there is significant scatter in the experimental 
data, the bulk of the observations are in agreement with 
the Hubbard-Douglas curve to within ^20%. The agree- 
ment with experimental data for circular cylinders falling 
horizontally and end-on (see Happel and Brenner 1965) is 
also good, to within ^10%. The cylinders falling end-on 
have a slightly higher fall speeds than the Hubbard & Dou- 
glas curve (2), whilst cylinders falling horizontally have 
a slightly lower fall speed than (2), which seems sensible 
since the Hubbard-Douglas theory is intended to represent 
the average over all orientations. 

Also shown in figure 3 is the MHKC prediction for 
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Figure 4. Normalised fall speed of randomly-oriented 6 arm bullet- 
rosette ice crystals (example shown inset). Black line shows the 
Hubbard-Douglas prediction as a function of bullet aspect ratio; 
grey line shows the corresponding MHKC prediction. Shading 
indicates fall velocities which lie outside the upper and lower 
bounds derived from analytical results (see text). 

randomly and horizontally oriented columns. As for pla- 
nar crystals, both MHKC curves overestimate the experi- 
mental data, by around 50% for an aspect ratio of 1.5, and 
by more than a factor of two for aspect ratios L/2a > 3. 
Upper and lower bounds may be constructed using equal- 
volume and enclosing spheres respectively as shown by 
the shaded areas. The MHKC calculations lie above the 
upper bound for all aspect ratios. 

2.3 Polycrystals 

Bullet-rosettes are often the dominant crystal habit in 
cirrus clouds (Heymsfield and Iaquinta 2000) which are 
typically composed of between two and six bullet crystals 
in a radial formation. Small 'embryonic' bullet-rosettes 
less than lOO^m in diameter have been observed (see 
figure 2 of Heymsfield and Iaquinta), and the fall speeds 
of these particles have been estimated using the Hubbard- 
Douglas and MHKC theories. There are, to the author's 
knowledge, no direct observations of the sedimentation 
velocity of these tiny rosettes. 

The capacitance of various bullet-rosette models has 
been calculated by Westbrook et al. (2008), who derived 
the fit C = 0.40/J x (L/2a)~ - 25 for rosettes with six 
arms (see inset figure 4). Here L is the length of the 
columnar section of the bullets, and (2a) is their width. 
The pyramid ends were assumed to be (L/2) in length. 
Using this data, the fall speeds of bullet-rosettes with a 
variety of aspect ratios has been calculated using equation 
2, as shown in figure 4. 

Also shown in figure 4 are the fall speeds of the 
same model crystals using the MHKC method. As for the 
capacitance calculations, random orientation is assumed. 
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At aspect ratios close to unity the MHKC fall speeds are 
~40% higher than the capacitance predictions, whilst for 
rosettes with thin arms L/2a = 4 the fall speeds are a 
factor of 3 larger than the Hubbard-Douglas prediction. 
This implies that the MHKC formula may be significantly 
overestimating the sedimentation rate of small rosette 
crystals, especially those with thin arms. This is confirmed 
by constructing an upper bound using an equal-volume 
sphere, as indicated by the shaded region. 

3 Discussion 

The sedimentation rate of small ice crystals with maxi- 
mum dimensions smaller than ~ lOO^rn have been esti- 
mated by exploiting the Hubbard-Douglas (1993) approx- 
imation for randomly-oriented particles of arbitrary shape, 
and the Roscoe (1949) approximation for planar particles 
settling horizontally. This appears to be the first time (to 
the author's knowledge) that these approximations have 
been applied to ice crystals. For hexagonal plate crystals 
the experimental data lies somewhere between the two 
approximations suggesting that the crystals were roughly 
horizontal in their orientation but with some flutter. For 
dendritic types the experimental scatter was much larger 
but appears to tell a broadly similar story. For columns 
and needles the sensitivity to crystal orientation seems to 
be weaker (given the fall velocities of circular cylinders 
falling horizontally vs. flat on) and the Hubbard-Douglas 
approximation gives a good estimate for the fall velocity 
of the crystals. The method has also been applied to model 
bullet-rosette crystals, for which there are no experimental 
data to compare with. 

The results are in substantial conflict with the MHKC 
parameterisation, which overestimates the fall speeds of 
these tiny crystals. From a physical point of view it is 
perhaps unsurprising that the MHKC boundary-layer for- 
mulation fails for Re< 1. However the suggestion in the 
past has been that equation 4 would hold down to these 
small sizes, largely based on the fact that the compar- 
ison with Stokes equation for a spherical particle is so 
good (see section 2). Unfortunately the comparisons with 
experimental data in this paper show that it does not hold 
in general for non-spherical ice crystals. Interestingly, 
Mitchell (1996) did make a comparison with the experi- 
ment of Yagi (1970) who observed ice crystals falling in 
a supercooled fog in Asahikawa, Japan. The average crys- 
tal fall speed was (v) = 0.107 ms -1 and the mean crys- 
tal size was close to 100^m. Mitchell used the recorded 
habit/diameter distribution, and mass-, area-diameter rela- 
tionships derived from other studies to estimate (v) = 
0.097ms~ 1 , in apparent agreement. However because the 
mass/area/size relationships had to be assumed from other 
studies, the comparison is rather uncertain; also the tiny 
fall speeds may be strongly affected by even the weakest 
up/downdraughts present within the fog. The comparisons 
made in this paper with the direct measurements of indi- 
vidual crystals by Jayaweera and Ryan (1972), Kajikawa 
(1973) and Michaeli (1977) in controlled laboratory con- 
ditions are a much stronger test of the MHKC formula, 
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and one which appears much less favourable to it at small 
Reynolds number. The upper bound calculations based 
on analytical solutions for insribed spheroids and equal- 
volume spheres confirm that equation 4 tends to provide 
an overestimate for the fall velocity of small crystals set- 
tling in a viscous air flow. 

It seems that capacitance is a good length scale to 
correlate ice particle fall speed data for very small crystals. 
However capacitance is not readily measured from in- 
situ probes (unlike A and D). For the crystal shapes 
considered here the capacitance may be calculated directly 
from D using the results in Westbrook et al. (2008). 
However there is evidence that many crystals in cirrus 
clouds may be complex polycrystals of irregular shape 
(Bailey and Hallett 2002, 2004). Further work is needed 
to estimate what capacitance such crystals are likely to 
have. The results given here however highlight the relative 
insensitivity to crystal shape, with almost all crystal types 
having a normalised velocity between 1 and 2. This 
indicates that even for rather tenuous, irregular crystal 
shapes (such as dendrites) the hydrodynamic radius R is 
not far removed from circular discs and spheres of the 
same overall dimensions. 

In this paper the focus has been on crystals whose 
dimensions have been measured explicitly in the lab, and 
we have found that the MHKC predictions generally over- 
estimated the measured fall velocity. Interestingly, A — D 
data from aircraft observations can nonetheless produce 
realistic fall velocities when applied to the MHKC param- 
eterisation in certain circumstances. For example Heyms- 
field and Miloshevich (2003) measured area ratios of 7 ~ 
0.75 ± 0.15 using replicator data collected in cold cirrus 
around — 60°C. This corresponds to normalised velocities 
in the range 1.3 to 2.0, and much of the experimental data 
is clustered toward the lower end of this range. The sug- 
gestion is that the crystals sampled were quite compact (as 
evidenced by the relatively high area ratio), a situation in 
which the MHKC formulation seems to perform best (eg. 
columns with L/2a ~ 1, figure 3). For more elongated or 
tenuous crystals typical of highly supersaturated or mixed- 
phase conditions, the evidence from figure 1, 2, 3 and 4 
is that the MHKC predictions are likely to be rather less 
realistic. 

For Re 3> 1 however, the MHKC relationship has 
proved very successful at accurately predicting ice parti- 
cle fall speeds, offering good agreement with experimen- 
tal data for most particle types. The experimental data 
presented here corresponds to Re < 0.4; it is not clear 
exactly where the cross-over from viscous to boundary- 
layer flow regimes occurs: in this work I have tentatively 
suggested a figure of ~ 100/mi based on the experimen- 
tal data. More work is needed to determine what happens 
at the intermediate Reynolds numbers which lie at the 
boundary between the two regimes. 

The fall speeds of the tiny crystals discussed in this 
paper are of interest for the interpretation of Doppler lidar 
measurements in cirrus. It may be possible to assess to 
what extent small crystals control the optical properties 
of cirrus cloud by examining the frequency distribution of 
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measured Doppler lidar velocities in such clouds. Using 
the results from section 2 we find that the maximum 
fall speed for a D = 60/mi crystal is 0.07ms -1 : this 
value corresponds to a compact hexagonal column with 
an aspect ratio of one, and a density equal to that of 
solid ice. The temperature was assumed to be — 60°C; at 
warmer temperatures the viscosity of the air is slightly 
higher and the crystals fall a little slower as a result. Note 
that the fall speed measurements by Jayaweera and Ryan 
(1972), Kajikawa (1973) and Michaeli (1977) included 
crystals over 100/im in size, and none of them fell faster 
than 0.06ms -1 . If the lidar backscatter is dominated by 
these tiny slow-sedimenting crystals, then the measured 
mean Doppler velocity should also be dominated by 
them since it is weighted by the backscatter. Analysis 
of 1 year of continuous measurements from a 1.5 /mi 
vertically-pointing Doppler lidar at the Chilbolton Facility 
for Atmospheric and Radio Research in Hampshire is in 
progress, and the results should help to inform the small 
ice debate. 

Finally the author notes that the Hubbard-Douglas 
approximation is likely to be applicable to non-spherical 
atmospheric aerosol particles, allowing estimates of the 
sedimentation rate and mobility of complex, irregular 
aerosols to be made. The only caveat is that the Hubbard 
and Douglas derivation assumes 'stick' boundary condi- 
tions, and so may not be suitable for sub-micron particles 
without an appropriate slip correction. 

4 Appendix A: randomly oriented crystals 

Here the average drag on a randomly oriented ice crystal 
is estimated using the approach of Hubbard and Douglas 
(1993). Consider the Oseen tensor, which describes a 
point hydrodynamic source (Happel and Brenner 1965): 

T < R ' = sb( I+ f) <« 

where I is the identity matrix, and R = r' r. The orien- 
tational average of T is given by: 

67777/2 ^ ^ 

where R = |R|, ie. T is identical to the Green's function 
for Laplace's equation (47ri?) _1 to within a constant 
factor. Based on this observation, consider the average flux 
of momentum away from the ice crystal: 

f aTds = -L x 0(r) (8) 
J 6nri 

where ds is a small element of surface area, and we 
define a to be the momentum flux density. Given the 
above considerations the scalar function 0(r) must satisfy 
Laplace's equation, and a may analogously be interpreted 
as the density of charge on the surface of a conductor with 
the same size and shape as the ice crystal. 
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The stress tensor S is constructed so as to conserve 
linear and angular momentum, and to ensure that u and 
the associated drag are co-linear: 

S = 6tt V [(V0)u o + u o (V0) - (V0) • (u )I] (9) 

where 4> = 1 on the surface of the particle, and cf> = in 
the far-field. Then the Stokes equations are: 

V ■ S = 6tt77U O V 2 = 0. (10) 

Essentially then, the flow system is described purely in 
terms of the electrical potential around a conductor of 
the same size and shape as the ice crystal. The angle- 
averaged fluid velocity at a given point is u(r) = uo^(r). 

Given the stress components above, the drag force on 
the particle is: 



67r?7 / S • nds 



(11) 



where n is the unit vector pointing normal to the crystal 
surface (outwards into the fluid). Using the analogy devel- 
oped above, the gradient may be interpreted in terms 
of the electrical potential gradient near a conductor of the 
same shape and size as the particle, ie (V0) sur f acc = — an. 
The total drag force is then simply: 



Drag = 67r?7|uo|C 
ie. Stokes formula (1) with R = C. 



(12) 



5 Appendix B: horizontally oriented planar crystals 

Here the drag on a thin (L <c 2a) planar crystal settling 
horizontally with velocity u is derived, following Roscoe 
(1949). The Stokes equations are: 



Y72 d P V72 

77V U — — — , TIN/ V 

OX 



— , 77V W 
dy 



dp 

dz 



(13) 



where the velocity components u,v,w (corresponding to 
the x, y, z axes) satisfy the continuity equation: 



du dv dw 
dx dy dz 



0. 



(14) 



and for convenience we assume u = v = w = on the 
crystal surface, and u — u ,v — 0,w — far from the 
crystal. The x direction is taken as the vertical axis. 
A flow system which satisfies (14) is: 



d(f> d(j) d(j) 

u = (p- x—, v = -x— , w = -x — 
ox dy dz 



(15) 



if V 2 = (ie. if cf> is a solution of Laplace's equation). 
Substituting this into equation 13 the velocities (15) satisy 
the Stokes equations provided that: 



P 



-277 



dx ' 



(16) 



Consider then an earthed conductor with the same shape 
as the crystal. Taking <j) = on the surface of the crystal 



and 4> — u Q far from the crystal matches the boundary 
conditions on u,v,w, and we identify (15) as the flow 
system round the ice crystal with cj) the electrical potential 
around the analogous conductor. 

The normal stress on the surface of the crystal is 
equal to p. The net force normal to the plate is therefore 
the difference between p on either side of the crystal, ie. 
Sp = 2?7 x 4-7rcr per unit area of surface (applying Gauss's 
law), where a is the charge density on the conductor. 
Integrating a over the whole surface leads to the total 
translational drag force on the crystal, ie: 



Drag = StttjuqC 



(17) 



since the total charge on the conductor is simply the capac- 
itance C multiplied by the applied voltage u . Comparing 
this result with Stokes law (1) leads to the conclusion: 



R=\ C 



(18) 



for planar ice crystals which are sufficiently thin relative 
to the dimensions of their horizontal cross-section. 
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